Absence of Dipole Glass Transition for Randomly Dilute Classical Ising Dipoles 
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Dilute dipolar systems in tliree dimensions are expected to undergo a spin glass transition as the 
temperature decreases. Contrary to this, we find from Wang-Landau Monte Carlo simulations that 
at low concentrations x, dipoles randomly placed on a cubic lattice with dipolar interactions do 
not undergo a phase transition. We find that in the thermodynamic limit the "glass" transition 
temperature Tg goes to zero as 1/y/N where A'^ is the number of dipoles. The entropy per particle 
at low temperatures is larger for lower concentrations {x = 4.5%) than for higher concentrations 
{x = 20%). 
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Disordered insulating materials often have randomly 
placed electric or magnetic dipoles that have long range 
dipolar interactions. Examples include impurities in al- 
kali halides that can be used for paraelectric cooling 
diluted ferroelectric materials jSi, disordered magnetic 
materials, and frozen ferrofluids y. These systems are 
typically modeled as spin glasses that have simpler in- 
teractions and yet are believed to capture the essential 
physics of interacting dipoles. Based on theoretical stud- 
ies of spin glasses with long range interactions [SIESIIIj 
one would expect dilute Ising dipolar systems to undergo 
a spin glass-like transition as the temperature decreases. 
However, in this paper we find the surprising result that, 
unlike the three dimensional Ising spin glass with 
interactions 0, dilute Ising dipolar systems do not 
undergo a phase transition as the temperature decreases. 
This may explain the lack of experimental evidence for 
such a transition in very dilute systems. 

An example of dipoles is two level systems (TLS) that 
dominate the physics of glasses at low temperatures ^3 ■ 
TLS often have randomly oriented electric dipole mo- 
ments that interact through an elastic strain field with 
a long range interaction that is a stress tensor general- 
ization of the vector dipolar interaction^^!- While there 
have been experimental hints of a spin glass transition 
among TLS in glasses at low temperatures 0| , there has 
been no definitive experimental proof that such a tran- 
sition occurs. Since the estimated concentration of TLS 
is low (100 ppm), our result may explain the absence of 
a transition even though TLS dipoles are randomly ori- 
ented and may not be Ising. 

Another example is the insulator LiHo2;Yi_j;F4 in 
which the holmium ions have Ising magnetic dipole mo- 
ments that lie along the z-axis due to crystal field effects 
Q. For very dilute systems (x = 4.5%) LiHoa;Yi_a:F4 
shows no sign of a transition [l5j . The lack of low tem- 
perature freezing in LiHoa;Yi_j;F4 has been attributed 
to dominant quantum mechanical effects in the so-called 
spin liquid or antiglass phase 0,^3- However, a theo- 
retical investigation of whether or not classical interact- 
ing dipoles undergo a spin glass phase transition at low 



concentrations has been lacking. Several previous stud- 
ies of dipolar interactions between randomly placed Ising 
dipoles have focussed on the ferromagnetic transition 
that occurs at higher dipole concentrations 0, 0, Il8j |. 
Monte Carlo simulations have looked at intermediate 
concentrations with x > 25% where there is a spin glass 
transition Xu et al. used mean field theory 

and found, depending on the lattice structure, ferromag- 
netic or antiferromagnetic transitions at higher concen- 
trations. They found a spin glass phase at lower spin 
concentrations but the properties of this phase were un- 
reliable because they had a replica symmetric solution. 
In short, there have been no definitive theoretical studies 
of the very dilute classical cases. In this paper we present 
the results of Wang-Landau Monte Carlo simulations on 
classical dilute Ising dipolar systems in three dimensions. 
We find that there is no phase transition for low concen- 
trations in qualitative agreement with experiment. 

In spin glasses the distribution P{q, T) of the over- 
lap order parameter q changes from being a Gaussian 
centered at g = at high temperatures to a bimodal 
distribution with peaks at g = ±1 at low temperatures. 
At intermediate temperatures it is relatively flat. We 
can define a characteristic glass transition temperature 
Tg as the temperature where P{q,T) is the flattest. In 
the thermodynamic limit we find that for a given dipole 
concentration Tg goes to zero as l/ViV where N is the 
number of dipoles. Also, we examine the entropy and find 
that for concentrations less than 20% there is a nonzero 
entropy per dipole as T — > 0. The entropy and lack of a 
transition are consistent with a large number of accessible 
low energy states and glassy behavior. 

We consider Ising dipoles randomly placed on a simple 
cubic lattice at concentrations of x = 4.5%, 12%, and 
20%. The interaction between any two dipoles pi and p2 
separated by a vector fi2 is given by the Hamiltonian: 



H{pi,P2) 



P1-P2- 3(fi2 • Pl)(fi2 • P2) 



(1) 



In addition to the energy units set by H, the units are 
set by Pi = ±z, the lattice constant a = 1, and Boltz- 
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mann's constant ks = I- These will be referred to as MC 
units where appropriate. We use the Ewald summation 
technique to handle the long range nature of the dipole 
interactions For x = 4.5%, the lattices had sites 
with L = 6, 8, 10, and 12, and for x = 12% and 20%, 
L — 4, 6, and 8. The number of dipoles N is the smallest 
even integer greater than or equal to xL^. 

The glassy energy landscape at low concentrations 
makes it difficult to equilibrate at low temperatures with 
the traditional Metropolis Monte Carlo approach. To 
overcome this, we have used the Wang-Landau (WL) 
Monte Carlo technique to calculate the density of 
states n{E) where E is the energy of the system. Briefly, 
this algorithm starts with an initial guess n{E) = 1 
and executes a weighted random walk on the energy 
landscape. Single flips of randomly selected dipoles are 
then accepted with a probability of min [1, n{Ei) / n(E f)] 
where Ei and Ef are the energies before and after the 
trial flip. If a step is accepted (rejected), then the den- 
sity of states is updated by the rule n{Ef(^i^) jn{Ef(^i)) 
where 7 > 1 is a scale factor. A histogram of the visited 
energies h{E) is recorded. The criterion for a satisfactory 
estimate of the density of states is given by the flatness 
of h{E), i.e., h{E) > e{h) for every energy E where 
< e < 1 determines the accuracy; typically, e « 0.95. 
Once the flatness condition is satisfied, the scale factor is 
set closer to 1 by the rule 7 — > ^/j, h{E) is reset to zero, 
and the algorithm is repeated. In all cases, we ran 20 
iterations with 7 starting at e and ending at 1.0000019, 
and n{E) was normalized such that X]_e'^(^) — 2^. 

The dipolar interaction is nearly continuous so each en- 
ergy bin may contain multiple states. We choose the bins 
to be as small as possible while maintaining reasonable 
computational times. The bin sizes depend on concentra- 
tion and system size. The bins are about 0.01 in units of 
energy per particle for 20% filhng and 0.001 for 4.5% and 
12%. The lowest temperature studied (T — 0.05) must 
be larger than the largest bin (0.02). We try to keep the 
bins small enough so that n{Eo) ~ 2, where Eq is the 
energy of the (degenerate) ground state. In all cases ex- 
cept one n{Eo) « 3.5. The exception (8'^ at 20%) has 
n{Eo) fa 10, so we discard the low temperature values of 
this system. 

We average over disorder by having different runs cor- 
respond to different quenched placements of dipoles with 
random initial orientations. The dipoles are fixed in po- 
sition but not in orientation. There are about 1000 runs 
for each x and L. As a check of our Wang-Landau pro- 
cedure, we were able to enumerate all the states for 1000 
different configurations for concentrations of 4.5% (L = 6 
and 8), 12% (L = 4) and 20% (L = 4) and determine the 
exact density of states. We found very good agreement 
with our WL results. 

Since we are looking for a spin glass phase, we define a 
generalized Edwards- Anderson overlap order parameter 
q = jj J2i Pi 'Pi 1 where is a dipole in the state of the 




FIG. 1: P{q) At X = 4.5%, L = 10 (46 dipoles). T =5, 1.6, 
1.1, 0.9, and 0.5. The lines transition from a Gaussian (T = 5) 
to a bimodal distribution {T = 0.5). 



current system, and p/ is a dipole in a low energy state 
found in a short, initial simulation |2^|2^. Then, to find 
the distribution P{q,E), q{E) is sampled and stored in 
a histogram during the simulation at the smallest scale 
factor where the estimate of the density of states is quite 
good. P{q, T) is calculated as 

P{q, T) = Ct{1/Z) n{E)P{q, E) exp{-E/kT) (2) 



where the sum is over all the energy bins and Z(T) = 
'^^n{E)eKp{—E/kT). Ct enforces normalization such 
that P{q, T) = 1 for every T. This method has been 
seen to give a reasonable order parameter distribution in 
the case of a Potts model |23| . 

It has often been convenient to find the spin glass 
transition temperature using Binder's g — [3 — 
(('?')/(9')')]/2, and {q^ = E,9"^^(9,T) 0. Since 
the ground state estimate is not the true ground state, 
we eliminate all runs in which g < 0.8 at the lowest tem- 
perature. If there is a second order phase transition, plots 
of g versus T for different size systems will cross at the 
transition temperature l^^. However, we find that these 
curves do not cross, so there is no second order spin glass 
phase transition. 

To investigate this further, we can look at how P{q, T) 
changes with temperature. For a system undergoing a 
phase transition, we expect P{q, T) to change from be- 
ing a Gaussian centered at 5 = at high temperatures 
to a bimodal distribution with peaks at g = ±1 at low 
temperatures. A typical example is shown in Figure ^ 
fora; = 4.5%,L = 10. We define a characteristic "glass" 
temperature Tg as the temperature where the distribu- 
tion P{q,T) is flattest. We define the deviation D{T) 
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Temperature (MC units) 

FIG. 2: Deviation from flatness of P{q, T) for x = 4.5% with 
sizes L = 6^, 8^, 10^, and 12^. The minima moves left with 
increasing size and defines a glass transition temperature. The 
plots for X — 12% and 20% are similar. 



from flatness in terms of the variance of P{q,T) as 

D{T) = L'(^[p{q,T)^(P{q',T))^,y'j (3) 

where {■ ■ ■)q indicates an average over all + 1 possi- 
ble values of q. D{T) is at a minimum when a plot of 
P{q,T) versus q is the flattest, defining Tg. A sample 
D{T) is plotted in figure El For a given dilute concentra- 
tion, Tg is tending to smaller temperatures as the system 
size increases which is consistent with — > as L — > cx). 
To find the dependence of Tg on the number N of dipoles 
for a fixed concentration, we plot the minimum of D(T) 
versus TV in figure 01 The best fit for dilute cases reveals 
that Tg ~ In contrast, D(T) for the ordered 

case (a; = 100%) yields a nonzero transition temperature 
independent of N. 

The absence of a transition is consistent with the ex- 
perimental finding that for very dilute systems {x = 
4.5%) LiH0a;Yi_3;F4 shows no sign of a transition [l5| . 
However, the absence of a transition in dilute dipolar sys- 
tems is unexpected since 3D Ising spin glasses with 1/r^ 
interactions undergo a phase transition P{q) for a 

spin glass and for a dilute dipolar system are different; in 
the thermodynamic limit as T ^ 0, P{q) for a spin glass 
has a few sharp peaks corresponding to ground state con- 
figurations separated by high barriers, while P{q) for the 
dilute dipolar system is fiat, indicating numerous accessi- 
ble low energy states separated by insignificant barriers. 
With very low barriers, states at both the top and bottom 
of the barrier contribute low energy states. The difference 
in barrier heights may be due to every site in a model spin 
glass being occupied so that in a spin glass with power 
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FIG. 3: Log-log plot of the maximum flatness for P{q) versus 
the number of dipoles at various concentrations. Open circles 
are 4.5%, squares are 12%, triangles are 20%, and upside down 
triangles are 100%. The left solid line has a slope of —1/2 
corresponding to Tg ~ N~^^'^. The right solid line has a 
slope of zero. Fits to Tg ~ iV"" yield a = 0.49(1) for 4.5%, 
a = 0.6(1) for 12%, a = 0.45(3) for 20%, and a = 0.02(5) for 
100%. The errors in the last digit are in parentheses. 



law interactions nearby spins will tend to have stronger 
interactions than distant spins and produce large barrier 
heights. In a dilute dipole system nearby sites are empty 
and so the low energy configurations are determined by 
distant dipoles which interact weakly and produce low 
barriers. 

The presence of many nearly degenerate accessible 
ground states is reflected in the flnite entropy per dipole 
near T = 0. We flnd that the low temperature entropy is 
larger for the lower concentration. We can calculate the 
total entropy Stot {T) directly from the density of states 
obtained by our WL Monte Carlo simulations: 

5,ot(r)-^^ + iogZ(r) (4) 

where {E{T)) is the average energy of the system. Stot is 
an absolute entropy and is not deflned relative to some 
reference value. To compare different system sizes, we 
consider the entropy per particle Sn = Stot/N where N 
is the number of dipoles. The entropy is very smooth, 
corresponding to a broad bump in the specific heat. 

To determine the entropy in the thermodynamic limit, 
we plot Sn{T) versus 1/N at a given temperature T. We 
fit a line to the data and then extrapolate to iV — > oo. 
Then, we plot the extrapolated value versus temperature 
(see Figure^. From Figure it is clear that the 4.5% 
and 12% cases have a nonzero entropy at low tempera- 
tures, but the 20% case is approaching zero. Finally, the 
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like phase transition as the temperature is lowered. This 
is consistent with having a large number of nearly de- 
generate accessible low energy states. Our result could 
explain the lack of experimental evidence for a transi- 
tion in LiHoa;Yi_2;F4 for small x and among two level 
systems in glasses at low temperatures. Thus, contrary 
to widely held notions, materials with dilute electric or 
magnetic dipoles cannot necessarily be modeled as spin 
glasses with long range interactions. 

We thank Stuart Trugman and Manoranjan Singh for 
helpful discussions. This work was supported by DOE 
grants DE-FG03-00ER45843 and DE-FG02-04ER46107. 

^Present address: The Salk Institute, La JoUa, CA 
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FIG. 4: Entropy per particle extrapolated to infinite size. 
From top to bottom the fillings are 4.5%, 12%, and 20%. The 
curve on the far right is for 100%. The error bars shown 
represent the standard error (~ 10"'^) of the distribution of 
S{T). Inset: Entropy at T = phase diagram. The hashed 
area is classically not accessible. The solid line is a guide to 
the eye. Above 20% ^(r = 0) is zero. 



extrapolated values are fit with a power law of the form 
AT^ + So, where A and A are constants, and So is a con- 
stant representing the zero temperature value of the en- 
tropy. The fit values are A = 1.1±0.2, A = 2.7±0.1, and 
So = (7.9±0.3) X 10-3 at 4.5%, A = 1.2±0.3, A = 2.9±0.2 
and So = (5.6 ±0.3) x 10"^ at 12%, and A = 0.37 ±0.05, 
A = 1.9 ± 0.1 and So = (-0.5 ± 0.5) x lO'^ at 20%. 
Note, the extrapolation at 20% gives a negative So, so 
it is zero; no actual data points have negative entropy. 
A phase diagram of the entropies at zero temperature 
is constructed in the inset of Figure ^ Notice that the 
low temperature entropy increases as the concentration 
decreases. This indicates that there are more accessible 
low energy states in systems with lower concentrations 
where the dipoles interact more weakly. Having a finite 
value of So implies that the zero temperature entropy 
5(0) may be nonzero, but this is not unprecedented for 
a classical system, e.g., nonintcracting spins. 

We do not think that the finite entropy near T = 
is due to the finite size of the energy bins. To test the 
effect of the bin size, we halved the bin size (doubled the 
number of bins) for the case of 6"^ at 20% and found an 
entropy change of about 5% which is consistent with the 
error estimates. We also ran the largest exact case (8'^ 
at 4.5%) through the WL algorithm with bins of width 
0.005, and found a change of 0.9% compared to the exact 
result with zero bin width. 

To summarize, we find the surprising result that at 
low concentrations {x < 20%) there is no spin glass- 
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